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ABSTRACT 


The solutions of the phase equation for a potential with the 
asymptotic properties of an atomic polarization potential are 
studied with the intention of developing a camputer program for 
direct integration of the phase equation for the s-wave phase 
shift of an arbitrary potential. Such a program is developed and 
its limitations are discussed. In addition, an appendix is devoted 
to the variational principles of Kohn and Schwinger for calculating 


phase shifts. 
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I. INTRODUCTION 


The Quantum Theory of scattering processes is of basic import- 
ance in studies of atomic, nuclear, and elementary particle physics. 
Analysis of Atomic collision processes is required in such diverse 
applications as gas dynamics, chemical reaction theory, and stellar 
dynamics. Nuclear Physics relies on scattering of various particles 
fram nuclei for a great deal of the quantitative information about 
nuclear processes while in Elementary Particle physics, very little 
other information is available except for the various scattering 
events that are observed. 

There are many approaches to the scattering problem, but one of 
the most useful (especially in low energy applications) is the method 
of partial wave analysis. The method was first introduced by Faxen 
and Holtzmark based on a method developed by Lord Rayleigh. Rather 
then delve into the basis of the theory, it will be assumed that the 
reader is familiar with the treatment of Schiff [Ref. 1] or Merzbacher 
(REE. 2). 

A particularly interesting procedure for calculating scattering 
phase shifts is the "Phase Method." Originally introduced by Courant 
and Hilbert [Ref. 3] and applied to scattering by Morse and Allis 
[Ref. 4], the method has recently been revived and expanded in a 
series of papers by Calogero [Ref. 5, 6] and his new book [Ref. 7]. 

This paper will investigate the application of the phase equation 
to a potential which has the asymptotic properties of an atomic 


polarization potential (1/x*) . The method employed is a direct 


integration of phase equation (11). In the appendix, the variational 
principles of Kohn and Schwinger are treated and are suggested as 


alternate methods for computation of phase shifts. 
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It... THEORY 


A. THE PHASE FUNCTION METHOD 

The phase function method is a method of reducing the second order 
linear Schrodinger Equation to a first order, non-linear Riccati equa- 
tion. For s-wave scattering, the solution for phase shifts reduces. 
to the solution of the differential equation 
(1) n(n) ee [very Je \ sant (xe +) 

To show how this equation. comes about, it is convenient to: approach 
the problem from the standpoint of A. Ronveaux [Ref. 8]. He starts by’ 
considering a second-order homogeneous differential equation with: two 


perturbative terms, Vy (X) and V. (X) ; 


(2) me + woul + Oida= Vdsu & KOO we’ 
with u(0) =p and u'(0) = v. This can be written in the form of the 


following system 


(3) Xie) = AGO Xx) + V(x) Xx) 
where 
ne On. Te 
X(x) (“| X(o) ‘? 
A(x) ?f O 1\ MVxjzf/O O 
Te -- ON, 


The unperturbed solution, xX, is, 


(4) X (x) = WON C 


aE 


where C = m 
C2 


is a constant vector, and 


uae a on Us 


ie 


and Uy and u, are two linearly independent solutions. The solution of 


the complete system may then be written 


(5) K(x) = Wx) C(x) 
Using the Lagrange Method of Variation of Constants, then 


OC, we } WOO) VO) Wix) CCx+} 


Cx) = M(x)C(x) 


CCo) = j wean fa) 


If we now define the phase function 


(6) Ses mac, Ge) / CK 


and differentiate, we find that 


Cn S Cadence + [Ma (x) 7 Hy (0) S SOx) = Mi(x)SO) 


where 


Evaluating M gives the phase equation 


S (x) - } J uc) | + | u, (x) so] } Jet VW/(x) 
x } Mia] us + Uy (x) S (x) | t V, (x) | u, (x) -ufcasealt 


where det(W)# 0 and 


{[-we-l Oy [+ [ucor]$ 


} | u ol u | + [ucry]f 


6) 
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Given the general form one now approaches the problem of s-wave 


phase shifts. The Schrodinger Equation (in natural units, oe | ) 


is 
(9) meee Ko eer) = Vor) adr 
choose 

Uo) ZT Wz O 

v.02). 


Then the general solution is 
Wir) 2 C, Sm ke ~- Co cos Kr 


~ Cae) seine — Sag) coskr | 
OF 
Wee Men) ces Wier Fy te sim (Ket) 


where 
S€q)= ¢, 13 4, ¢r) = ~- tan 4, (rpx) 
Then, letting V, = Bio eG Chup sinks, U, = cosk ‘ga clef W: K the general 


phase equation becomes 
7a 


(io) =(d/dr) tan Yo == [ very /k] (sim kr + tan yocoskr | 
with tan), (0,k)=0, Noting that (d/dr)tann, = sec ye(dyo/dr} 


it becomes apparent that 


(11) Cd/dr) No = I vor) /k | sin’ (Kr + "0 ) 


B. THE POTENTIAL FUNCTION 

The object of this paper is to investigate the properties of the 
numerical solutions for the Phase Equation (11) when applied to a 
potential which has the asymptotic properties at infinity of an 


Atomic Polarization potential, i.e., 


Ce) a (const. ) x ifr d 
[<7 oe 


Je 


In order to simplify the analysis, an additional condition was imposed 
on the potential in the form of a requirement that as [-? 0 , the 


potential should approach zero as ( , l.e., 


Clee 





> (const xr 
T——9 © 


The function chosen to satisfy these conditions was 


¥ r—27 oo 4 
o 


The asymptotic properties of this form are obvious, however, it is 


f ay 
UO ae a dee 6), ae eek 


difficult to calculate with. In order to simplify the problem further, 


the substitution pe xe was chosen which leaves 
{ 


| : i rx 
(12) Cine a F lorena 
o 
It is no longer readily apparent that U(r) has the correct asymptotic 


properties, but it may be easily verified by integration and power series 


expansion. The success of this form is shown in Fig. l. 
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Figure 1 


Plot of the Potential Function 
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III. CALCULATIONS 


Having derived the phase equation (11), one now wishes to devise 
a method to put it to practical use. Calogero develops a variational 
scheme for the phase equation [Ref. 9] while other authors deal with 
the equation in other ways [Ref. 10, 11]. An excellent review has 
been prepared by Babikov [Ref. 12] which provides a good general 
introduction to actual calculations and an equally good bibliography. 
Owing to the fact that the equation is non-linear, the various ap- 
proaches have dealt with same form of approximate solution. This 
paper will deal with the relatively simple problem of constructing 
a program to directly integrate the phase equation using numerical 
approximation only. 

The particular method chosen was a fourth order Runge-Kutta 
Integration scheme to be implemented on the IBM-360 time sharing 
system. General discussions of the numerical method may be found 
in numerous places [Ref. 13, 14], and so only the necessary portion 
of that theory is presented here. 

The Runge-Kutta method proposes to solve the equation 
i) v a f ( x, Vv ) 
subject to the initial conditions X = Xo when Y = Yo° The iterative 
solution process proceeds via the following equation; 

(14) vi - vy, + nel 4+ 2k, + 2KRaee Ky | 
© 
Kas re eee) 
Ke? h {Cx nt 2 
Ky2 hf 0X4 2h) Wt a &) 


where 
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Nig Mf Cknt hy Yn + Ks) 


and where X,,,, = X,¢h . The programming needed to provide a 
solution of the equation (11) which gives the exact phase shift, 

uy BMis ee Ny ©; k) , then reduces to programming a simple iterative 
procedure. The program developed is printed at the back and is in 
such a form that it may be run either on the time-sharing system for 
quick approximations or directly on the computer system in the regular 
batch processing. A second version of the program used to obtain a 
print out of the values of the phase function at various r values is 
also included. The use of the second form will become obvious as 

the discussion unfolds. 

It now becomes necessary to consider the actual running of the 
program and its particular limitations. In order to properly investi- 
gate the program, the Yukawa potential, 

~mr 
VOr)# 4 c / mr 
is convenient to consider since Calogero treats the case with strength 
parameter, g = -10, m=1. The output produced by the second program 
is graphed in Figs. 2 and 3. 

The most obvious feature of the graphs is the presence of "levels" 
which arise from the presence of possible bound states in the potential. 
Calogero discusses the connection of this with Levinson's theorem and 
shows that this structure in the phase function can be used to locate 
bound states. This very feature though poses a distinct problem to 
the evaluation of the exact phase shift. Owing to the presence of this 
structure, a simple test for the end of the integration that depends on 


the phase function tending to a constant value at infinity is useless. 
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Figure 2 


Plot of Phase Function for K = 0.1, Yukawa Potential 
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Plot of Phase Function for K = 1.0, Yukawa Potential 
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For this reason, the integration must be carried out past any possible 
structure before a cutoff of the above type is allowed. In order that 
the program will do just that, an integer counter is established and 
required to be larger than a certain lower limit which is designed 
to carry the integration over the structure. The procedure for doing 
this is as follows; the second form of the program is first run with 
the intention of finding the structure and determining the value of DR 
which gives the best results. Having found the structure, the counter 
may be set appropriately. By use of this artifice, the program is 
enabled to reproduce the results given in Calogero's book to within 
a few percent. It should be noted that another way of treating this 
problem is to establish a lower limit on the basis of the value of r 
reached. This approach would be superior if there were a method for 
predicting the position of the last flat region in No" The counter 
method was chosen only on the basis of a desire to know the number 
of iterations used, for computer timing purposes. 
A further pecularity with the direct integration routine is the 
sensitivity to the s MRED, DR. In investigations with different 
40 ra values of DR, it was found that when 
DR is chosen too large the structure 
at small r tends to destroy the ac- 
20 x curacy of the routine by adding 
¥ large terms over some of the flat 
x regions in the phase function. An 


| 3 
pA example is for K = 0.1 and dr = 0.1. 


Figure 4 The result is es 44.1 in 110 steps 


as shown in the diagram (Fig. 4). 
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The true result (in agreement with Calogero's book, p. 22) is 

De 6.24 in 38 steps found by taking dr = 0.01 or less. Now the 
utility of the second program for preliminary tests become apparent 
since it allows us to establish the number of bound states for the 
potential and the proper increment and the lower counter limit in 
just a few test runs which may be quickly run on the computer time 
sharing system. 

The program also has a very important limitation in that the 
accuracy decreases sharply for a small k due to numerical round off 
errors arising from the necessity of dividing by a small number. The 
Situation may be improved by choosing a smaller increment, DR, but 
after a certain point this procedure tends to become self defeating 
because decreasing DR increases the program running time. A test run 
on a Yukawa potential with g = -2, K = 0.1, and dr = 0.1 gives a close 
answer to the zero energy phase shift, but if the value of k is reduced 
to K = ion) the accuracy in the calculated value of the phase shift 


actually decreases. 


Ja 


IV. RESULTS 


The final problem is the application of the programs to the 
polarization potential discussed earlier, eqn. (12). Figures 5 
and 6 detail the results of the explicit calculation with the 
second program. They indicate that a value of DR of 0.1 will carry 
the program nicely through to the asymptotic region if the lower 
counter limit is set at 200. The lack of structure in the graphs 
tends to make this requirement useless though and it becames pos- 
Sible to set the limit at 10 iterations just to let the program get 
started before the asymptotic cutoff procedure is allowed to function. 
The lack of structure just commented on tends to indicate a lack of 
bound states in the potential as might be suspected from the very 
weak nature of the function. Tests could be done by imposing a 


strength parameter on the potential, i.e., 
{ 


Vir) = qe |dxx%e 


O 
which could then be adjusted as suited the investigator. By means 


t X 


of such an artifice, bound states might begin to appear when the peak 
of the potential is high enough. In addition, increasing the strength 
of the potential would increase the magnitude of the phase shift. One 
might expect that such structure would appear for values of r less than 
1.0 and that the way to locate such structure would be to decrease the 
value of the increment DR. Though interesting, such an investigation 
does not have any apparent physical use since the only physical feature 
of the potential is the asymptotic dependence at large r. Thus further 


discussion of this point has been amitted here. 
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Figure 5 


Plot of Phase Function for K = 0.1, Polarization Potential 


23 


Figure 6 


Plot of Phase Function for K = 1.0, Polarization Potential 
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The most interesting feature of the calculation fram the stand- 
point of possible application of the program to real, physical 
potentials is the necessity of carrying out the calculation to larger 
values of the radial coordinate than was the case for the Yukawa 
potential. The reason for this behavior is obvious, the polarization 
potential only drops off as fast as Ve" at infinity. This leads to 
the observed difference in the convergence of the scheme and probably 
to the reason why the only places that the Phase method is applied in 
the literature are in papers dealing with nuclear type potentials which 
may be assumed to fall off very rapidly. To get some idea of the magni- 
tude of this effect, the graph in Fig. 2 may be examined in comparison 
to the graph in Fig. 5. The plot for the Yukawan potential has sta- 
bilized by r = 1 while with the polarization potential, the plot is 
just stabilizing at r = 25. These observations hold great import for 
any application of the method to atomic like potentials since the matter 
of time available for computation is often very critical. Even on the 
IBM-360 time sharing system, timing for the program had become critical 
when increments less than 0.1 were chosen. In addition, output disk 
space became a problem. The problem of timing could be expected to be 
considerable worse on smaller computing systems. 

Having discussed the phase function and its properties, use may 
now be made of program 1 to obtain a plot of ny VS- energy as amieeeyateye. 7 « 
The most apparent feature of Fig. 7 is the 0.78 radian peak at k = 0.36. 
It is a real peak as opposed to an apparent peak that could be produced 
by a poor choice of the increment at such a low value of k. An inter- 
pretation of this peak will not be attempted at this point but a com 


ment is appropriate on the behavior of Ny (kK) on the low energy portion 


Z> 


of the peak. It seems to drop off drastically, possibly to zero. 
Below k = 0.1, no data taking was attempted due the necessity of 


integrating too far with too small an increment. 
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Figure 7 


Plot of Phase Shift vs. K for Polarization Potential 
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V. CONCLUSIONS 


The program developed for the phase equation integration has been 
tried on a test case and applied to a second, semi-physical case with 
somewhat the expected results. This program has in these two cases 
performed well as long as the precautions mentioned were taken. It 
should be noted that no error analysis of the routine was attempted 
due to lack of time, however, the comparison of the results on the 
Yukawa potential for g = -10 with the results published by Calogero 
in his book [Ref. 7], tends to point to an accuracy that is within 
a few percent. Whether this result is valid for any potential other 
than the Yukawa potential and in particular whether it is equally 
true for atomic level potentials is a question that remains unan- 
swered at this point. It is expected that simple numerical errors 
arising from the Runge-Kutta integration scheme and the ordinary 
computer errors of round off, etc. will decrease the accuracy of the 


calculation. 
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APPENDIX 


(A.1) The Kohn-Hulthen Method 
The most frequently used variational procedure is the Kohn- 
Hulthen method [Ref. 15, 16] which is based on the stationary property 
of the functional ~ \ where 
(Al) se am |YCE-H) bdr 
Le 
and A is same function of the phase shift determined fram the normal- 


ization of b as —-700° . This paper will employ tangent normalization 





throughout, i.e., 
(A2) p 7 Ga ° | je (Kr) + tang ne Cer) | 
r—> co 


where Ne is the phase shift for the ee 


partial wave and fe and vy 
are the spherical Bessel and Neuman functions respectively. Consider- 


ing the special case of s-wave scattering, Q = 0, then (A2) becames 
Y a, (da)? [ sinke/ir A. Pay Yo MM: /cx) 
J ; 
< (uz) * | _ sin (Kirt 7) 


COSY, RG 





With this form for p , ohne now considers the explicit variation of (Al) 
«(AT jsucv 1K Ud AT + vests ce. U)Sedr 
Since for s-wave, p has no angular dependence, 

be BEC H)Y 


so that the second integral in (A3) becomes 


\+\ 3 ("Ged + Kr - v| od ar 


raday\ dr 
- fan _— aba" Ate 4 \6 dco) ode 
yr 


Zo 


Integrating by parts, this gives 


“an [as fd {oe age} - 4 _ dd *d (AY) 


dr 


1 (shee -vo) de 


da \or} Weed (Sv) - 5¥r? at) 8d (gy 


4 a VAT 


é ts (a a rac oe srs | be jshenwita 
, at Jr wn 


SI - 2 (sboviexs aybar 4 Ja[e(vash. adds 


The last term, upon application of the boundary conditions, reduces 


a 5 Yan yo /k } so that 
sy (I+ tan Ye /k) BZ 150s =) ee 


Thus (A3) becomes | aa 


But the Schrodinger equation is just 
(Vi4e-v)P=0 


so that 


(A5) dCi + ton No /k ) : 


The stationary functional, [A] = Tan Yio /K , is given by 


(A6) (oe We ee a 2 \d(e- wd bdr 
An 
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The trials function, p , 1s chosen to have the form 
(A7) W = a” 4 lige X; (0) 


where uy = (di) “aes K i | and the X; are same convenient basis 
functions, that rapidly approach zero as r-”?@ . The asymptotic part 
of the wave function a , 1s chosen so as to shield the = singularity 
in cos kr/r as r~7 0. It is evident that p has the correct asymptotic 


form as 


| 
Y = (ui \? [sinkr/er ~ \k Csinkr /p2p? - coskr /er)| 


aa 


> (uz) 


ns [— 





(sin br /ker) + | (cos Kr/r) | 


7 eee 


The Kohn-Hulthen principle then states that the variation of {\] with 


respect to all arbitrary constants, l.e., 
(a8) SI\l = S(A+I)=- 0 


produces a value of [\] that is second order accurate ind) el he 


explicit variation is done by writing the functional in the form; 


0 A] 


i 


at; City + 226, @ + B+ 
zy a) 
where 


2 


MM; = 2m \ x. CE- WX; a Mii 


\ Xj (e-wiP 


A 
sP 


B= 2m \ oe-m¢ 
AS 


Sil 


Variation of the n constants, C, then gives the system of linear equations 
(A10) z Mi Ci = = Ry 
J 


Since the trial function is linear in \ A 
(0) C 


%; = ic; + \ I? ; 


and since the /1, | are independent of \ 


COF Cc! 
Co ce + \ C; 


This leaves the equations: 
Co) 
2. Mi le = ie; 
Ci) 
= IZ; 


(All) 


ve 


2 My P C/ 
Solving for the ir,.8 : 


TM = wy + Aw, + dw, 4 Bo + \8, 4) Bat) 


where 
a Hiss (wo) Co) (0) 
W, : Vs oF ‘ M;; C +t 2 2 C; Ie, 
e | d 
Cr) = Ci) =e 
— - ¢ Co) : 
WW = py Ch S45 C. a + 2 2 G Ky 16 Z é C, Kk, 
3) é 


Ci iy yy Ci) € 
26° My Gf? + 2 26fPR 
J yA 


Then a final variation with respect to \ gives, to first order 


\- = Guinea. ae 
2S; + Gs) 


and the second order accurate, stationary functional is 


PAS ae ech [hy = fw + B +1 a Bit!) 
ee SU, + By) 
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(A.2) The Schwinger Method 

A less frequently used but more elegant variational method is 
that due to Schwinger [Ref. 17, 18, 19, 20). It gives the stationary 
value of the functional [Kcoty | directly. The problem is to solve 
the radial Schrodinger Equation. 


(A13) Au + a = Bc2a) | ug() = Or) ur) 
oh i“ 


To do so, one solves for the Green's Function Giger it ) given by, 


ene i 
which gives 


me) gh Garin + [at aenaan | Galen Stree 


Ge (rir, )= ~K iro paler) my kre) GS la 
(A1L5) 
Se O(c Cig emai aire 


The integral equation for ?(r) is then 


a> 
(Al16) Din lee ) = ae ep = fdr. GGr tn, ) OG arene) 


For large r, the asymptotic forms are 
CIC ORS (rr- 2 (0+1)7) 


> sin(yr - 4 (es1)q) 





Kr yoCicr ) 





K¢ Ne (Kr) 


therefore, the asymptotic form of JU,(r) may be written as 





? yous ee err 
O> 


+ sin eps : oa] | de Vy fe) uCry / 
© 
But, this is to match to the asymptotic form 


U pg ( op 
(A17) 


UgCr) —-> Cos Ikr - “(e+1)i | + sinlkr- 3 (esi ] tanne 
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and so one sets 


(eos ee] 


SPR “{e i \ae. te yelX to) OCr,) Ug 6 ‘a0 


Inserting the integral form for yecer) in terms of Ue(r) fram (Al6) 


ae (i> 7 
IS 


therefore, 
Ktanv\g 5 \3: DARE: \ p Vbey Or) ugte) GOI) Ul) Ug) 
© © Oo 
dividing now by to ns "Ne 


oa 


ur) 


\1 Uor)ug dr) + {o |e Ur ue ce y Gale, JOOD Uely) 


ev 


(Als) K cot 2 5 


Bos A “cee Oe re ee te ee 8 ee me ee a ow ee Oe 


en 2 
| \a- if car yunerdjecer | 


.@) 


For s-waves this becames 


Vale GUE je fr Or duce) Grin) Ulry) olre) 
(ALS) Keedyy = vie | ee yee 


Te ee : 
i \sr Vic er yeaa cr | 


where G(rir,) is now: 


| 
Cigars ) ea aay Cos Kr, Ge Oe 
= u Sn KT pcos le res £ 
K 


In order to show that (Al19) is in fact stationary to second order for 


the exact wave function, consider 


(A20) S\ kot v0 | = | keed 7 (ot du)| = [kcot (uf 
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Expanding and collecting terms, we have 
S(kcot no | > (ar UCr)o(r) | eed - (ar Gr V¥y) Uy) ote | 


Gt 


+ 2 Ar Gy Wow fucry- (cert. ye WL. 2 


| le dy Uiervirsntr | +2] [i [ar oer Sosinte [xf Jor ver UCrolr) si mat | 
6 lO 5 


yD Oo© 
(a. Ur) vcr) Juco - |r. Gr try) Uerruens | 
gO eines ° 


: me 
\: dy Ul (r)olr) sin ky | 
Cry 


combining the terms the numerator for the expression will be 


i 2 So CK? at 
|. (ar UL sors | Nt UG Ge Jen jas Gr UGH ¥0%)| 
“ ), & v 


+ 2 jar Ola duce) Je = js Fy Gr | rye) Or, ) UCf,) \\ 
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Equation. 
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solution for the C$. This result is shown in Kohn's Paper (Ref. 15]. 
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